purpose is to highlight an approach to spatial expression analysis using ideas from interactive visualization
high-level idea: both interactive and direct visualization give us something like maps. We want to be able to inspect them simultaneously.
We want to generate these interactive views in demand, so we wrap the core visualization code in an htmlwidgets package (https://www.htmlwidgets.org/).
First, we"ll install some packages needed to prepare the data and run the analysis.
We use spatial data analysis packages to prepare the geojson which ends up as input to our visualization
library("dplyr")
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library("geojsonio")
##
## Attaching package: 'geojsonio'
## The following object is masked from 'package:base':
##
## pretty
library("raster")
## Loading required package: sp
##
## Attaching package: 'raster'
## The following object is masked from 'package:dplyr':
##
## select
## Registered S3 method overwritten by 'geojsonlint':
## method from
## print.location dplyr
data_dir <- file.path(Sys.getenv("HOME"), "Data") download_data(data_dir)
## Warning in dir.create(directory, recursive = TRUE): '/github/home/Data' already
## exists
## [[1]]
## [1] "/github/home/Data/mibiSCE.rda"
##
## [[2]]
## [1] "/github/home/Data/masstagSCE.rda"
##
## [[3]]
## [1] "/github/home/Data/TNBC_shareCellData"
loaded_ <- load_mibi(data_dir) cell_full <- read_csv(file.path(data_dir, "TNBC_shareCellData", "cellData.csv"))
## Parsed with column specification:
## cols(
## .default = col_double()
## )
## See spec(...) for full column specifications.
input_data <- function(cell_full, tiff_paths, sample_id, crop_size=1000, simplify=0.4) { # read the associated raster im <- raster(tiff_paths[grepl(paste0("p", sample_id, "_"), tiff_paths)]) im <- crop(im, extent(im, 0, crop_size, 0, crop_size)) # comment out to work on full data polys <- polygonize(im) %>% filter(!(cellLabelInImage %in% c(0, 1))) # filter down to relevant data cell_data <- cell_full %>% filter( SampleID == sample_id, cellLabelInImage %in% polys$cellLabelInImage ) polys <- polys %>% left_join(cell_data) # run U-Map on channel data channels <- setdiff(colnames(cell_data), c("cellLabelInImage", "SampleID", "cellSize", "Background", "tumorYN", "cell_cluster", "tumorCluster", "Group", "immuneCluster", "immuneGroup")) dimred <- umap(cell_data[, channels]) %>% .[["layout"]] %>% data.frame() %>% rename(V1 = X1, V2 = X2) dimred <- cbind(cell_data, dimred) # convert to geojson geo_polys <- geojson_json(polys) %>% geojson_sp() %>% ms_simplify(keep = simplify) %>% geojson_json() list("dimred" = dimred, "polys" = geo_polys) }
The tumor clusters (from the original rda object) are <span style="background:#d9d9d9">4</span>, <span style="background:#bc80bd">7</span>, <span style="background:#ccebc5">10</span>, <span style="background:#ffed6f">17</span>.
The immune clusters are <span style="background:#8dd3c7">1</span>, <span style="background:#ffffb3">2</span>, <span style="background:#bebada">3</span>, <span style="background:#fb8072">4</span>, <span style="background:#80b1d3">8</span>, <span style="background:#fdb462">10</span>, <span style="background:#b3de69">11</span>, <span style="background:#fccde5">12</span>.
data_ <- list() for (i in seq_len(6)) { data_[[i]] <- input_data(cell_full, loaded_$tiffs, i) } linked_views(data_[[1]]$polys, data_[[1]]$dimred)
linked_views(data_[[2]]$polys, data_[[2]]$dimred)
linked_views(data_[[3]]$polys, data_[[3]]$dimred)
linked_views(data_[[4]]$polys, data_[[4]]$dimred)
linked_views(data_[[5]]$polys, data_[[5]]$dimred)
linked_views(data_[[6]]$polys, data_[[6]]$dimred)